##########################################################################################

library(data.table)
library(optparse)
library(ArchR)

##########################################################################################
option_list <- list(
    make_option(c("--comine_data_file"), type = "character"),
    make_option(c("--out_path"), type = "character") 
)

if(1!=1){
    
    ## 整合atac和rna的文件
    comine_data_file <- "~/20231121_singleMuti/results/subcell/cluster17/cluster17.combineRNA.motif_peak2gene.Rdata"

    ## 输出
    out_path <- "~/20231121_singleMuti/results/motif_position"

}

###########################################################################################
parseobj <- OptionParser(option_list=option_list, usage = "usage: Rscript %prog [options]")
opt <- parse_args(parseobj)
print(opt)

comine_data_file <- opt$comine_data_file
out_path <- opt$out_path

dir.create(out_path , recursive = T)

###########################################################################################
## 导入数据
a <- load(comine_data_file)
## atac_proj

###########################################################################################
## 细胞类型
cell_type <- unique(atac_proj@cellColData$cell_type)
cluster <- unique(atac_proj@cellColData$seurat_clusters)

if( length(which(cluster == c(12,13))) == 2 ){
    cluster <- 17
}

###########################################################################################
# Get locations of motifs of interest:
motifPositions <- getPositions(atac_proj, name="Motif")
motifGR <- stack(motifPositions, index.var="motifName")
motifGR <- data.frame(motifGR)
motifGR$cell_type <- cell_type
motifGR$cluster <- cluster

out_file <- paste0(out_path, "/cluster" , cluster , ".motif.Rdata")
save( file = out_file , motifGR )